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Abstract 

The eiectronic structure near a singie classical magnetic impurity in a super- 
conductor is determined using a fully self-consistent Koster-Slater algorithm. 
Localized excited states are found within the energy gap which are half elec- 
tron and half hole. Within a jellium model we find the new result that the 
spatial structure of the positive- frequency (electron-like) spectral weight (or 
local density of states), can differ strongly from that of the negative frequency 
(hole-like) spectral weight. The effect of the impurity on the continuum states 
above the energy gap is calculated with good spectral resolution for the first 
time. This is also the first three-dimensional self-consistent calculation for a 
strong magnetic impurity potential. 
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Magnetic impurities have a dramatic effect on superconductivity. Most work on the 
experimental and theoretical effects of magnetic impurities have focused on bulk thermody- 
namic quantities, such as the reduction of the superconducting transition temperature T c 
with increasing magnetic impurity concentration Theoretical approaches to treat the 
thermodynamic effects of impurities include Born scattering (Abrikosov-Gor'kov theory) |2j 
and approximate solutions to all orders in the impurity potential for classical [|||-[5|] and quan- 
tum || spins. A key issue addressed in Ref. 3 was the evolution of localized excited states 
around impurities into an impurity band, eventually leading to gapless superconductivity at 
high enough impurity concentrations. The effects of impurities on bulk properties have also 
been treated within a strong-coupling formalism, but not self-consistently or beyond the 
Born approximation (e.g. Ref. 7). Concern about bulk properties in the above and related 
work did not extend to properties very near the impurity. 

Among the first properties calculated in the vicinity of an impurity were the structures 
of screening clouds around a charged impurity and a magnetic impurity PJTOfl in a 
superconductor (characterized by exponentially-decaying Friedel-like oscillations). The os- 
cillation of the order parameter around a magnetic impurity was first evaluated without 
self-consistency [|TT|,|i"2"|] . A self-consistent calculation of the order parameter at the impurity 



and very far away for weak impurity potentials was done by Schlottmann [|T^]. Interest in lo- 



cal properties near impurities has been revived by advances in scanning tunneling microscopy 



(STM) near impurities embedded in a metallic medium ||T4 . 

Motivated by the possibility of measuring the local electronic structure near an impurity 
in a superconductor with an STM, the differential conductivity through an STM tip was 
calculated near impurities for superconductors with isotropic and anisotropic gaps within the 



Born approximation . The differential conductivity measured at a point x and voltage V 
and temperature T can be related to the local density of states (LDOS) at the tip location 
as follows: 

dI(x,V,T) 1 , dn(uj) /lmG(x,x;w = eV)\ 

dv = w L duj ^r{ n )■ (1) 



Here e is the charge of the electron, n{u) is the Fermi occupation function at temperature 
T, N Q is the density of states at the Fermi level, and dl/dV is in units of the normal-metal's 
differential conductivity. The local density of states is the imaginary part of the retarded 
Green's function fully dressed by the interaction of the electronic system with the impurity. 
Within the Born approximation, then, the differential conductivity can be expressed in 
terms of the retarded homogeneous Green's functions of the superconductor, g(x, x'; uj) and 
/(x, x'jtu) (where / is the anomalous Green's function), in the following way: 
dl(x, V) 



oc Im 



ff(x,x; eV) 

(2) 



dV 

+<?(x, 0; eV)V g(0, x; eV) - /(x, 0; eV)V f(0, x; eV) 



where the impurity is at the location in the solid and the non-magnetic impurity potential 
is Vq5(x). For a BCS superconductor with an isotropic order parameter in a parabolic band, 
for uj much smaller than the Fermi energy, 

<?(x, x'; u) = -^e-^^^r/^A ( cog fc + u sin kpr ) 
k F r \ VA 2 - uj 1 ) 

/(x, x'; uj) = ^Le-«^ A sin k F r (3) 

where r = |x — x'|, A is the order parameter, and ^ = kkp/^A is the coherence length. 
These expressions are valid for uj above and below A so long as the imaginary parts of both 
/ and g are multiplied by sgncu. 

The Born approximation calculation would not yield any localized states around the 
impurity, and does not consider the effect of the change in electronic structure on the local 
order parameter A(x). Recent preliminary STM results indicate certain features of the local 
density of states near a magnetic impurity fllEf : (1) discrete states are evident within the 
energy gap and the LDOS associated with them is asymmetric with voltage around V — 0, 
and (2) the LDOS becomes indistinguishable from the bulk density of states within a distance 
greater than the Fermi wavelength but much less than the coherence length. A calculation 
beyond the Born approximation where the order parameter is self-consistently determined 
and the continuum spectrum is found would be useful for explaining these results. 



We will present a fully self-consistent calculation of the local electronic structure near a 
magnetic impurity which is based on a Koster-Slater-like Green's function technique rather 
than the Bogoliubov-de-Gennes (BdG) equations [17]]. The BdG equations are Schrodinger- 
like equations for the electron and hole components of the quasiparticle wavefunctions. The 
localized states of a vortex core in a superconductor and other properties revealed by the 
LDOS were calculated self-consistently using the BdG equations |18|J19[| . These results were 
compared with a measurement of a single vortex by an STM on superconducting NbSe2 
|2(|. The BdG equations were successfully used again to explain the STM observations of 



the vortex lattice |2TJ . 

Since the original application of the Green's function algorithm to localized vibrational 



modes |[22|| , it has been applied to numerous problems including conduction electrons in 
metals [23[] and impurity states in magnets ||24|| , but not to superconductors. To place this 
algorithm in context we will contrast it with the BdG equations. To find the quasiparticle 
wavefunctions, typically the inhomogeneity is placed in a sphere of radius R with appropriate 
boundary conditions. The value of R is determined by the spectral resolution necessary 
for accurately determining the order parameter self-consistently and the spectral width of 
features measurable by (for example) the STM. The typical complications resulting from 
approximating an infinite system by a finite-size system apply, such as discrete states above 
the energy gap and the heavy investment of computer time required for large values of R. 

The Green's function method works within a sphere whose radius is determined by the 
range of the inhomogeneous potential. In essence we shall invert the Dyson equation in 
real space. The Dyson equation for an inhomogeneity in a superconductor, the Gor'kov 



equations [25], can be written as: 



J dx" [5(x - x") - g CT (x, x", w)V(x")] G CT (x", x; to) = g CT (x, x; 



(4) 



where 
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the same relationship exists among Go-, G a and F a , and 



V(x") 



( aV s {x") + V {x") <5A(x") ^ 



(6) 



^ 5A(x") <rV s (x") - Vb(x") ^ 
5A(x) = A(x) — A D , where A(x) is the position-dependent order parameter and A Q is the 
order parameter of the homogeneous superconductor. Vs(x) is a localized, spin- dependent 
potential, such as one originating from a classical spin. Vo(x) is a localized non-magnetic 
potential. The quantization direction of the superconductor's spins (cr = ±1/2) is chosen 
parallel to the direction of the classical spin. For the purposes of heuristic arguments in 
this Letter we will assume ferromagnetic coupling between impurity spin and electron spin 
(Vs < 0). Antiferromagnetic coupling only produces the trivial change that the antiparallel 
spin quasiparticle is attracted to the classical spin. g a is given by Eqs. (|3]) with N Q labelling 
the density of states for each spin. G a ^ G- a due to the differences in the potential in 
Eq. (H). The self-consistency equation is 

/oo 
do;n(a;)ImF (T (x, x; uo) — A Q , (7) 

where 7 is the pairing strength. 

One strength of this formalism is its reliance on the short-range nature of the inhomo- 
geneous potential. Solution of Eq. (£|) requires inverting the frequency-dependent real-space 
matrix M(u;) = <5(x — x') — g CT (x, x'; cj)V(x'). For x"s where V(x') is negligible, that portion 
of M(u;) is triangular and trivially invertible. Hence the numerical complexity of inverting 
M(c<j) is governed by the range of V(x). The radius of the sphere the system is solved in is 
governed by the range of the longest-range potential, which in this paper will be 5A(x). 

As a special case of Eq. (|), the energies of the localized states within the gap correspond 
to those u > where detM(u;) = 0. We model the impurity with a Gaussian Ks(x) of range 
kp x . Fi gure 1 shows the dependence of the localized state energies for the first two angular 
momentum channels on the strength of the magnetic potential, v s = 7tN \a J dxVs(x)|. For 
Vq = 0, the Eqs. (§)-(0) are unchanged by the transformation u> —u and s — > — s so 



the poles of G a (u) must come in symmetric pairs (cr = ±1/2) around to — p6 |. The 



quasiparticle state for small v s corresponding to these poles consists of an electron (u > 0) 
in the spin band parallel to the classical spin, which we will label up (f ), and a hole (u < 0) 
in the other spin band (J,). A hole in the down band has spin up and is attracted to the 
classical spin. Spin is a good quantum number for the quasiparticle, which is spin up. The 
spatially-integrated spectral weight is one-half each for the electron and hole components 
(this is true for all localized quasiparticle states in this model). 

For small v s an already existing quasiparticle can be bound by the classical spin. This 
allows for local pair-breaking excitations of energy less than 2A D . As the potential strength 
increases, the excitation energy of each angular momentum state £ decreases, at some critical 
value (v* e ) vanishes, and then begins to increase. For v s > v* the ground state contains a spin- 
up quasiparticle bound to the classical spin ||. The low-energy excitation now corresponds 
to removing that quasiparticle, which is a spin-down excitation. This qualitative behavior 
of the excitation energies can be extracted from the Shiba model [|[, where the magnetic 
potential is modeled by a delta function at the impurity site. In Ref. 3 zeros of M(uj) are 
found, neglecting the component of Reg(x, x'; uj) which is symmetric for uj — > —uj. These 
poles are shown in Figure 1. 

Even in the normal state the spectral weight of the up band has a peak at the origin, while 
that of the down band is pushed away from the origin. In the superconductor this asymmetry 
is evident in the quasiparticle spectral weights. Figure 2 shows the calculated angular 
momentum I = spectral weights in the up and down bands for two values of v s , indicating 
that the asymmetry of the spatial structure of the spectral weights becomes more pronounced 
with increasing v s . Also shown are I = 1 spectral weights for v s = 0.875. The asymmetric 
localized-state spectral weights produce differential conductivities which are asymmetric as 
well |[27| . Figure 3 shows differential conductivities for v s = 1.75 at three locations — right 



at the impurity, kp 1 away, and lOkp 1 away. By lOkp 1 the spectrum has recovered to the 
homogeneous spectrum. The differential conductivity farther away can be easily recovered by 
constructing the self-consistent T-matrix V(x')<5(x' — x") + V(x')G(x', x"; uj) V(x") directly 
from the dressed Green's functions. The spatial dependence of the differential conductivity 



far away from the impurity is qualitatively identical to the Born result ]T5 . 



The asymmetric behavior of the electron-like and hole-like spectral weights due to the 
difference between the up band states and down band states does not emerge from the model 
of Ref. 3, which implies identical spatial dependences of the spectral weights for the up and 



down bands. As pointed out by Koster and Slater [23|, proper treatment of the symmetric 
real part of g{0, 0; uj) is essential for obtaining the local electronic properties around an 
impurity The proper approximation for the real symmetric part of the Green's function is 
to average it over a small volume given by the range of the potential modeled by the delta 
function. We therefore introduce a new parameter into the Shiba model, a, which is the 
a;-symmetric part of Re(g(0, 0; uj = 0))/ttN o . The brackets indicate that g{0, x; uj = 0) has 
been averaged over a small volume. For a^O the ratio of the spectral weight of the up band 
at the impurity, Af(0), to the spectral weight of the down band at the impurity, A^{0), is 

A T (0) 1 + 2av s + (1 + a 2 )v 2 s 



i 2 \ 2 



(8) 



A^O) 1 - 2av s + (1 + a 2 )v 

The introduction of a does not change the localized state energies qualitatively. 

We now discuss the structure of A(x). Figure 4 shows <5A(x)'s for two values of v s . 
While oscillating with wavelength ~ nkp 1 , <5A(x) falls off to a negligible potential within 
lOkp 1 . This justifies performing the numerical calculation in a sphere of that radius. A 
typical radial grid of 100 points provides a numerically robust solution. The self-consistent 
solution also depends on the value of N A /kp, which for niobium is 3.6 x 10~ 5 p8 



As shown in Figure 4, for large values of v s , A(x = 0) < 0. Sign changes in A, as seen 
in pair tunneling, have been suggested for magnetic impurities in the barriers of Josephson 



junctions [29|. Our sign change in A(0) occurs (at T = 0) precisely at v* . The symmetry 
of Eqs. (f|)-([7p under u) — > —uj and s — > — s implies that ImF^r, r,w) = — ImFj (r, r, — uj) . 
As the pole in the spin-up band goes from electron-like (u > 0) to hole-like {uj < 0) and 
the pole in the spin-down band goes from hole-like to electron-like the contribution to A(0) 
changes sign abruptly. A(0) as a function of v s is shown in the insert of Figure 4. At T > 
the transition would be smoothed somewhat. 



The behavior of A(0) as a function of v s comes from the introduction at v* of a quasi- 
particle into the ground state of the system. The spin up quasiparticle localized near the 
impurity in the ground state suppresses the local order parameter. Exciting the low-energy 
state (removing the spin-up quasiparticle) causes a pair to enter the condensate. We find 
that exciting the low-energy state for v s > v* increases A (0), whereas excitation of quasipar- 
ticles typically reduces A(x) (which is the case for v s < v* . Also, exciting the low-energy 
state reduces the induced spin of the superconductor at the impurity. We note that this 
picture should have implications for the theory of gapless superconductivity, which is now 
based on the formation of impurity bands through hybridization of single-quasiparticle ex- 
cited states around impurities, without consideration of the coherence of excitations with 
the condensate. 
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FIGURES 

FIG. 1. Solid lines indicate the self-consistently calculated localized excited state energies for 
angular momentum channels 1 = and 1. The analytic model of Ref. 3 is shown by the dashed 
lines. At a critical value of v s = v* , the character of the I = state changes from spin up to spin 
down. The kink evident in the solid line is real, and is due to the discontinuous change (at T = 0) 
in A(x) at v* . 

FIG. 2. Spectral weights for the i = state's up-band component (electron-like for v s < v*) 
and down-band component (hole-like for v s < v*) as a function of position for two values of the 
magnetic potential. Also shown are the spectral weights for the up and down bands for i = 1 for 
v s = 0.875. 

FIG. 3. Differential conductivity (dl/dV), relative to the normal metal at three distances (in 
units of kp 1 ) from the impurity for v s = 1.75, showing the evolution of the spectrum from one 
dominated by the localized states near r = to one dominated by the bulk spectrum at r = 10. 
The curve at r = has been shrunk by a factor of two so that it can appear on the same scale 
as the other two curves. The presence of the large peak near the origin on the negative- frequency 
side of the spectrum indicates that the impurity is strong (v s > v* ). dl/dV has been evaluated 
from Eq. (1) with [3 = 7.5/A G , which for niobium corresponds to T ~ 2K. 

FIG. 4. Change in the local order parameter, 5A(x), for two values of the magnetic potential. 
The change becomes negligible beyond lOA;^ 1 . (Insert) Order parameter at the impurity, A(0), as 
a function of v s , indicating a discontinuous change at v* . 
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